Bienvenid@s a la primera tarea del curso Statistical Thinking. Esta tarea tiene como objetivo evaluar los contenidos teóricos de la primera parte del curso, los cuales se enfocan principalmente en análisis exploratorio de datos y conceptos introductorios de probabilidades. Si aún no han visto las clases, se recomienda visitar los enlaces de las referencias.
La tarea consta de ejercicios prácticos con el fín de introducirlos a la programación en R enfocada en el análisis estadístico de datos.
Slides de las clases:
Videos de las clases:
En la siguiente sección deberá resolver cada uno de los experimentos computacionales a través de la programación en R. Para esto se le aconseja que cree funciones en R, ya que le facilitará la ejecución de gran parte de lo solicitado.
Para esta pregunta usted deberá trabajar en base al conjunto de datos
hearth_database.csv, el cual esta compuesto por las
siguientes variables:
En base al dataset propuesto realice un análisis exploratorio de los datos (EDA). Para su análisis enfoquen el desarrollo en las siguientes tareas:
Respuesta
my.frame <- read.table(file = "hearth_database.csv",header = T,sep = ",", fileEncoding = "latin1")
colnames(my.frame)[colnames(my.frame) == "ï..age"] <- "age"
#slope
mean(my.frame$slope)
## [1] 1.39934
median(my.frame$slope)
## [1] 1
quantile(my.frame$slope)
## 0% 25% 50% 75% 100%
## 0 1 1 2 2
max(my.frame$slope)
## [1] 2
#ca
mean(my.frame$ca)
## [1] 0.7293729
median(my.frame$ca)
## [1] 0
quantile(my.frame$ca)
## 0% 25% 50% 75% 100%
## 0 0 0 1 4
max(my.frame$ca)
## [1] 4
#thal
mean(my.frame$thal)
## [1] 2.313531
median(my.frame$thal)
## [1] 2
quantile(my.frame$thal)
## 0% 25% 50% 75% 100%
## 0 2 2 3 3
max(my.frame$thal)
## [1] 3
#age
mean(my.frame$age)
## [1] 54.36634
median(my.frame$age)
## [1] 55
quantile(my.frame$age)
## 0% 25% 50% 75% 100%
## 29.0 47.5 55.0 61.0 77.0
max(my.frame$age)
## [1] 77
#trestbps
mean(my.frame$trestbps)
## [1] 131.6238
median(my.frame$trestbps)
## [1] 130
quantile(my.frame$trestbps)
## 0% 25% 50% 75% 100%
## 94 120 130 140 200
max(my.frame$trestbps)
## [1] 200
#chol
mean(my.frame$chol)
## [1] 246.264
median(my.frame$chol)
## [1] 240
quantile(my.frame$chol)
## 0% 25% 50% 75% 100%
## 126.0 211.0 240.0 274.5 564.0
max(my.frame$chol)
## [1] 564
#thalach
mean(my.frame$thalach)
## [1] 149.6469
median(my.frame$thalach)
## [1] 153
quantile(my.frame$thalach)
## 0% 25% 50% 75% 100%
## 71.0 133.5 153.0 166.0 202.0
max(my.frame$thalach)
## [1] 202
#oldpeak
mean(my.frame$oldpeak)
## [1] 1.039604
median(my.frame$oldpeak)
## [1] 0.8
quantile(my.frame$oldpeak)
## 0% 25% 50% 75% 100%
## 0.0 0.0 0.8 1.6 6.2
max(my.frame$oldpeak)
## [1] 6.2
#summary(my.frame)
cor(my.frame[,7:14])
## slope ca thal age trestbps
## slope 1.00000000 -0.08015521 -0.10476379 -0.16881424 -0.12147458
## ca -0.08015521 1.00000000 0.15183213 0.27632624 0.10138899
## thal -0.10476379 0.15183213 1.00000000 0.06800138 0.06220989
## age -0.16881424 0.27632624 0.06800138 1.00000000 0.27935091
## trestbps -0.12147458 0.10138899 0.06220989 0.27935091 1.00000000
## chol -0.00403777 0.07051093 0.09880299 0.21367796 0.12317421
## thalach 0.38678441 -0.21317693 -0.09643913 -0.39852194 -0.04669773
## oldpeak -0.57753682 0.22268232 0.21024413 0.21001257 0.19321647
## chol thalach oldpeak
## slope -0.004037770 0.386784410 -0.57753682
## ca 0.070510925 -0.213176928 0.22268232
## thal 0.098802993 -0.096439132 0.21024413
## age 0.213677957 -0.398521938 0.21001257
## trestbps 0.123174207 -0.046697728 0.19321647
## chol 1.000000000 -0.009939839 0.05395192
## thalach -0.009939839 1.000000000 -0.34418695
## oldpeak 0.053951920 -0.344186948 1.00000000
#boxplot(my.frame[,7:14],main="Hearth")
boxplot(my.frame$slope,main="Hearth slope")
boxplot(my.frame$ca,main="Hearth ca")
boxplot(my.frame$thal,main="Hearth thal")
boxplot(my.frame$age,main="Hearth age")
boxplot(my.frame$trestbps,main="Hearth trestbps")
boxplot(my.frame$chol,main="Hearth chol")
boxplot(my.frame$thalach,main="Hearth thalach")
boxplot(my.frame$oldpeak,main="Hearth oldpeak")
hist(my.frame$age[my.frame$target == "YES"],main="Distribución de Edad para 'YES'")
hist(my.frame$age[my.frame$target == "NO"],main="Distribución de Edad para 'NO'")
hist(my.frame$slope[my.frame$target == "YES"],main="Distribución de slope para 'YES'")
hist(my.frame$slope[my.frame$target == "NO"],main="Distribución de slope para 'NO'")
hist(my.frame$ca[my.frame$target == "YES"],main="Distribución de ca para 'YES'")
hist(my.frame$ca[my.frame$target == "NO"],main="Distribución de ca para 'NO'")
hist(my.frame$thal[my.frame$target == "YES"],main="Distribución de thal para 'YES'")
hist(my.frame$thal[my.frame$target == "NO"],main="Distribución de thal para 'NO'")
hist(my.frame$trestbps[my.frame$target == "YES"],main="Distribución de trestbps para 'YES'")
hist(my.frame$trestbps[my.frame$target == "NO"],main="Distribución de trestbps para 'NO'")
hist(my.frame$chol[my.frame$target == "YES"],main="Distribución de chol para 'YES'")
hist(my.frame$chol[my.frame$target == "NO"],main="Distribución de chol para 'NO'")
hist(my.frame$thalach[my.frame$target == "YES"],main="Distribución de thalach para 'YES'")
hist(my.frame$thalach[my.frame$target == "NO"],main="Distribución de thalach para 'NO'")
hist(my.frame$oldpeak[my.frame$target == "YES"],main="Distribución de oldpeak para 'YES'")
hist(my.frame$oldpeak[my.frame$target == "NO"],main="Distribución de oldpeak para 'NO'")
Pruebe el teorema central del limite aplicando un muestreo de la media en las distribuciones Gamma, Normal y una a su elección. Grafique los resultados obtenidos y señale aproximadamente el numero de muestreos necesarios para obtener el resultado esperado, pruebe esto con las siguientes cantidades de muestreo \(\{10,100,1000,5000\}\). ¿El efecto ocurre con el mismo número de muestreo para todas las distribuciones?.
Respuesta
# Definición de variables o estructuras necesarias para el muestreo.
# Primero definimos la función plot_gauss, que nos permitirá graficar rápidamente.
plot_gauss <- function(x, name) {
avg <- mean(x)
std <- sd(x)
x_axis <- seq(min(x), max(x), length = length(x))
y_axis <- dnorm(x_axis, mean = avg, sd = std)
hist(x, main = sprintf("Gráfico de la distribución %s", name))
lines(x_axis, y_axis, type = "l", col = "red")
}
# Representamos las distribuciones Normal, Gamma y de Poisson, de parámetros N(1.5, 0.75), G(2, 1.5) y P(5).
avg <- 1.5
std <- 0.75
shape <- 2
scale <- 1.5
lambda <- 5
# Inicializamos los vectores que acumularán las medias para cada muestreo.
iter <- 1000
normal_vector <- vector(length = iter)
gamma_vector <- vector(length = iter)
poisson_vector <- vector(length = iter)
# Graficamos las 3 distribuciones para distintos tamaños de la muestra
sample_size <- 10
# Realizar el muestreo de las distribuciones.
for(i in 1:iter) {
new_vector <- rnorm(n = sample_size, mean = avg, sd = std)
normal_vector[i] <- mean(new_vector)
new_vector <- rgamma(n = sample_size, shape = shape, scale = scale)
gamma_vector[i] <- mean(new_vector)
new_vector <- rpois(n = sample_size, lambda = lambda)
poisson_vector[i] <- mean(new_vector)
}
# Código para gráficos
par(mfrow = c(1, 3))
plot_gauss(normal_vector, sprintf("Normal con tamaño %s", sample_size))
plot_gauss(gamma_vector, sprintf("Gamma con tamaño %s", sample_size))
plot_gauss(poisson_vector, sprintf("de Poisson con tamaño %s", sample_size))
# Graficamos las 3 distribuciones para distintos tamaños de la muestra
sample_size <- 100
# Realizar el muestreo de las distribuciones.
for(i in 1:iter) {
new_vector <- rnorm(n = sample_size, mean = avg, sd = std)
normal_vector[i] <- mean(new_vector)
new_vector <- rgamma(n = sample_size, shape = shape, scale = scale)
gamma_vector[i] <- mean(new_vector)
new_vector <- rpois(n = sample_size, lambda = lambda)
poisson_vector[i] <- mean(new_vector)
}
# Código para gráficos
par(mfrow = c(1, 3))
plot_gauss(normal_vector, sprintf("Normal con tamaño %s", sample_size))
plot_gauss(gamma_vector, sprintf("Gamma con tamaño %s", sample_size))
plot_gauss(poisson_vector, sprintf("de Poisson con tamaño %s", sample_size))
# Graficamos las 3 distribuciones para distintos tamaños de la muestra
sample_size <- 1000
# Realizar el muestreo de las distribuciones.
for(i in 1:iter) {
new_vector <- rnorm(n = sample_size, mean = avg, sd = std)
normal_vector[i] <- mean(new_vector)
new_vector <- rgamma(n = sample_size, shape = shape, scale = scale)
gamma_vector[i] <- mean(new_vector)
new_vector <- rpois(n = sample_size, lambda = lambda)
poisson_vector[i] <- mean(new_vector)
}
# Código para gráficos
par(mfrow = c(1, 3))
plot_gauss(normal_vector, sprintf("Normal con tamaño %s", sample_size))
plot_gauss(gamma_vector, sprintf("Gamma con tamaño %s", sample_size))
plot_gauss(poisson_vector, sprintf("de Poisson con tamaño %s", sample_size))
# Graficamos las 3 distribuciones para distintos tamaños de la muestra
sample_size <- 5000
# Realizar el muestreo de las distribuciones.
for(i in 1:iter) {
new_vector <- rnorm(n = sample_size, mean = avg, sd = std)
normal_vector[i] <- mean(new_vector)
new_vector <- rgamma(n = sample_size, shape = shape, scale = scale)
gamma_vector[i] <- mean(new_vector)
new_vector <- rpois(n = sample_size, lambda = lambda)
poisson_vector[i] <- mean(new_vector)
}
# Código para gráficos
par(mfrow = c(1, 3))
plot_gauss(normal_vector, sprintf("Normal con tamaño %s", sample_size))
plot_gauss(gamma_vector, sprintf("Gamma con tamaño %s", sample_size))
plot_gauss(poisson_vector, sprintf("de Poisson con tamaño %s", sample_size))
Realice el experimento de lanzar una moneda cargada 1000 veces y observe el comportamiento que tiene la probabilidad de salir cara. Para realizar el experimento considere que la moneda tiene una probabilidad de \(5/8\) de salir cara. Grafique el experimento para las secuencias de intentos que van desde 1 a 1000, señalando el valor en que converge la probabilidad de salir cara.
Respuesta
# Inicializa un vector para almacenar las probabilidades en cada paso
probabilidades <- numeric(1000)
# Simular lanzamientos
for (lanzamientos in 1:1000) {
# Genera una secuencia de lanzamientos de la moneda y cuenta cuántas caras salieron
resultados <- sample(c("cara", "sello"), lanzamientos, replace = TRUE, prob = c(5/8, 1 - 5/8))
caras <- sum(resultados == "cara")
# Calcula la probabilidad de salir cara en este punto
probabilidad_actual <- caras / lanzamientos
probabilidades[lanzamientos] <- probabilidad_actual
}
# Gráfico de la convergencia
plot(1:1000, probabilidades, type = "l", col = "blue", xlab = "Número de intentos", ylab = "Probabilidad de salir cara")
abline(h = 5/8, col = "red", lty = 2)
legend("topright", legend = c("Probabilidad teórica", "Probabilidad simulada"), col = c("red", "blue"), lty = 2:1)
Remontándonos en la televisión del año 1963, en USA existía un programa de concursos donde los participantes debían escoger entre 3 puertas para ganar un premio soñado. El problema del concurso era que solo detrás de 1 puerta estaba el premio mayor, mientras que detrás de las otras dos habían cabras como “premio”.
Una de las particularidades de este concurso, es que cuando el participante escogía una puerta, el animador del show abría una de las puertas que no fue escogida por el participante (Obviamente la puerta abierta por el animador no contenía el premio). Tras abrir la puerta, el animador consultaba al participante si su elección era definitiva, o si deseaba cambiar la puerta escogida por la otra puerta cerrada. Un vídeo que puede ayudar a comprender mejor este problema en el siguiente link.
Imagine que usted es participante del concurso y desea calcular la probabilidad de ganar el gran premio si cambia de puerta en el momento que el animador se lo ofrece. Utilizando listas/arrays/vectores simule las puertas del concurso, dejando aleatoriamente el premio en alguna posición del array. Hecho esto, genere un numero de forma aleatoria para escoger una de las puerta (posiciones de la estructura), para luego ver si cambiando de posición tendrá mayores posibilidades de ganar el premio. Genere N veces el experimento y grafique cada una de las iteraciones, tal como se hizo en el ejercicio de las monedas.
Respuesta:
# Creamos una función que simule el juego
montyhall <- function(cambiar = TRUE){
puertas <- sample(1:3,3) #Puertas donde la posición que tiene el 3 es el premio
eleccion <- sample(1:3,1) #Elección del participante.
puertas_no_elegidas <- setdiff(puertas, eleccion)
puertas_no_elegidas_corr <-puertas_no_elegidas[puertas_no_elegidas != 3]
if(length(puertas_no_elegidas_corr) == 1){
puerta_abierta <- puertas_no_elegidas[1]
}
else{
puerta_abierta <- sample(puertas_no_elegidas[puertas_no_elegidas != 3], size = 1)
}
eleccion_final <- eleccion == 3
if(cambiar){
puertas_restantes <- setdiff(puertas, c(eleccion, puerta_abierta))
eleccion_final <- puertas_restantes
eleccion_final <- puertas[eleccion_final] == 3
}
return(eleccion_final) # Retornamos la elección, esta puede que tenga el premio o no
}
# Función que simula N juegos
n_juegos <- function(n = 10 ,cambiar_puerta = TRUE){
resultados <- replicate(n, montyhall(cambiar_puerta))
probabilidad_cambio <- mean(resultados)
plot(1:n, cumsum(resultados)/(1:n), type = "l", col = "blue",
xlab = "Número de simulaciones", ylab = "Probabilidad de ganar ")
abline(h = probabilidad_cambio, col = "red", lty = 2)
legend("topright", legend = c("Probabilidad", "Probabilidad promedio"), col = c("red", "blue"), lty = 2:1)
}
n_juegos(1000)
n_juegos(1000,FALSE)
Ustedes disponen de los dados D1 y D2, los cuales no están cargados y son utilizados para comprobar que \(\mathbb{P}(AB)=\mathbb{P}(A)\mathbb{P}(B)\) cuando el evento A es independiente del B. Para estudiar la independencia considere que los eventos A y B se definen de la siguiente manera; sea A el evento dado por los valores obtenidos en el lanzamiento del dado D1, este está compuesto por \(A=\{D1=1,D1=2,D1=6\}\). Por otro lado, el evento B viene dado por los valores obtenidos con el dado D2, el que está conformado por \(B=\{D2=1,D2=2,D2=3,D2=4\}\). Con esto, tendremos un \(\mathbb{P}(A)=1/2\), \(\mathbb{P}(𝐵)=2/3\) y \(\mathbb{P}(AB)=1/3\). Compruebe de forma gráfica que al realizar 1000 lanzamientos (u otro valor grande que usted desea probar) se visualiza que \(\mathbb{P}(AB)=\mathbb{P}(A)\mathbb{P}(B)\).
Hecho lo anterior, compruebe el comportamiento de un segundo grupo de eventos, dados por el lanzamiento de solo el dado D1. Donde, los eventos para D1 quedan definidos como: \(A =\{D1=1,D1=2,D1=6\}\) y \(B=\{D1=1,D1=2,D1=3\}\). ¿Se observa independencia en este experimento?.
Se le aconseja que para simular los lanzamientos de dados utilice la
función sample() para generar valores aleatorios entre 1 y
6. Compruebe los números generados por la función con los casos
favorables de cada uno de los eventos a ser estudiados.
Respuesta:
n_lan <- 6000000
# Primer grupo de eventos
dice1 <- sample(1:6, n_lan, replace = TRUE) # Dado 1
dice2 <- sample(1:6, n_lan, replace = TRUE) # Dado 2
dice_count1 <- table(dice1)
dice_count2 <- table(dice2)
favorable_a <- c(1, 2, 6) # {1, 2, 6}
favorable_b <- c(1, 2, 3, 4) # {1, 2, 3, 4}
favorable_ayb <- c(1, 2) # {1, 2}
count_a <- sum(dice_count1[favorable_a]) # |A|
count_b <- sum(dice_count2[favorable_b]) # |B|
# Ahora contamos los casos donde A y B se cumplen al mismo tiempo
count_ayb <- 0 # |A∩B|
for (i in 1:n_lan) {
if (dice1[i] %in% favorable_a && dice2[i] %in% favorable_b) {
count_ayb <- count_ayb + 1
}
}
# Calculamos las probabilidades
prob_a <- count_a / n_lan # P(A)
prob_b <- count_b / n_lan # P(B)
prob_ayb <- count_ayb / n_lan # P(A∩B)
print(sprintf("P(A) = %s", prob_a))
## [1] "P(A) = 0.4999925"
print(sprintf("P(B) = %s", prob_b))
## [1] "P(B) = 0.666550833333333"
print(sprintf("P(A)*P(B) = %s", prob_a*prob_b))
## [1] "P(A)*P(B) = 0.333270417535417"
print(sprintf("P(A∩B) = %s", prob_ayb))
## [1] "P(A∩B) = 0.333211833333333"
Ahora visualizamos las probabilidades
par(mfrow = c(1, 2))
y_ax <- c(0, 0, 1, 1)
x_ax <- c(0, 1, 1, 0)
# Graficamos P(A)
plot(1, ann = FALSE, type = "n", axes = FALSE, xlim = c(0, 1), ylim = c(0, 1), main = "P(A)")
x_p_a <- c(0, prob_a, prob_a, 0)
x_p_a_c <- c(prob_a, 1, 1, prob_a)
polygon(x_p_a, y_ax, col = "blue")
text(mean(x_p_a), mean(y_ax), expression(P(A)), col = "white", cex = 1.5)
polygon(x_p_a_c, y_ax, col = "red")
text(mean(x_p_a_c), mean(y_ax), expression(P(A^c)), col = "white", cex = 1.5)
# Graficamos P(B)
plot(2, ann = FALSE, type = "n", axes = FALSE, xlim = c(0, 1), ylim = c(0, 1), main = "P(B)")
y_p_b <- c(prob_b, prob_b, 1, 1)
y_p_b_c <- c(0, 0, prob_b, prob_b)
polygon(x_ax, y_p_b, col = "blue")
text(mean(x_ax), mean(y_p_b), expression(P(B)), col = "white", cex = 1.5)
polygon(x_ax, y_p_b_c, col = "red")
text(mean(x_ax), mean(y_p_b_c), expression(P(B^c)), col = "white", cex = 1.5)
Ahora grafiquemos \(\mathbb{P}(AB)\)
plot(1, ann = FALSE, type = "n", axes = FALSE, xlim = c(0, 1), ylim = c(0, 1), main = "P(A∩B)")
polygon(c(0, prob_a, prob_a, 0), c(prob_b, prob_b, 1, 1), col = "blue")
text(prob_a / 2, 1 - prob_b / 4, expression(P(A * "∩" * B)), col = "white", cex = 1.5)
polygon(c(0, prob_a, prob_a, 0), c(0, 0, prob_b, prob_b), col = "purple")
text(prob_a / 2, prob_b / 2, expression(P(A * "∩" * B^c)), col = "white", cex = 1.5)
polygon(c(prob_a, 1, 1, prob_a), c(prob_b, prob_b, 1, 1), col = "purple")
text(1 - prob_a / 2, 1 - prob_b / 4, expression(P(A^c * "∩" * B)), col = "white", cex = 1.5)
polygon(c(prob_a, 1, 1, prob_a), c(0, 0, prob_b, prob_b), col = "red")
text(1 - prob_a / 2, prob_b / 2, expression(P(A^c * "∩" * B^c)), col = "white", cex = 1.5)
Notemos que el área del cuadrado de la esquina superior izquierda, que representa el área de \(\mathbb{P}(A \cap B)\), tiene la misma área que la superposición de las figuras anteriores,
Un amigo ludópata suyo le comenta que el truco de jugar en el casino esta en no parar de apostar y apostando lo mínimo posible. Ya que así, tienes mas probabilidades de ganar el gran pozo que acumula el juego. Usted sabiendo la condición de su amigo, decide no creer en su conjetura y decide probar esto a través de un experimento.
Para realizar el experimento usted decide asumir las siguientes declaraciones, bajo sus observaciones:
En el primer experimento deberá obtener la evolución de los fondos hasta que el jugador se queda sin fondos para jugar. Puede ser útil seguir la lógica de una moneda cargada para realizar esto. Pruebe esto con una apuesta igual a 5, 25 y 50 graficando los resultados. Comente los resultados obtenidos.
Para la segunda parte de este experimento, con las funciones ya creadas, proyecte 5000 apuestas y obtenga la probabilidad a la que converge el experimento empezando con una apuesta de: 5, 25 y 50. Para este experimento considerar como éxito (1) cuando su función ruina supera los 200 y considere lo contrario cuando su función perdida es menor o igual a 0.
Respuesta
# Función para obtener el desarrollo de las apuestas
ruina <- function(money = 100, bet = 5) {
vec_money <- c(money)
while (0 < money && money < 200) {
game <- as.integer(runif(1, 1, 19))
sign <- if (game <= 8) 1 else -1
money <- money + sign * bet
vec_money <- append(vec_money, money)
}
return(vec_money)
}
plot(ruina(), type="l", col="blue", xlab="N° de juegos", ylab="Fondos", main="Evolución de los fondos (apuesta = 5)")
plot(ruina(bet = 25), type="l", col="blue", xlab="N° de juegos", ylab="Fondos", main="Evolución de los fondos (apuesta = 25)")
plot(ruina(bet = 50), type="l", col="blue", xlab="N° de juegos", ylab="Fondos", main="Evolución de los fondos (apuesta = 50)")
En general se ve que la tendencia es perder, solo que al apostar el mínimo se pierde más lento ya que al ser una variación más pequeña, se desestabiliza (es decir, se va a los extremos) de a poco, no así al apostar 50, que es posible perder en dos juegos.
n_bet <- 5000
bet5_vector <- vector(length = n_bet)
bet20_vector <- vector(length = n_bet)
bet50_vector <- vector(length = n_bet)
for (i in 1:n_bet) {
final5_money <- tail(ruina(bet = 5), 1)
final20_money <- tail(ruina(bet = 20), 1)
final50_money <- tail(ruina(bet = 50), 1)
win5 <- if (final5_money <= 0) 0 else 1
win20 <- if (final20_money <= 0) 0 else 1
win50 <- if (final50_money <= 0) 0 else 1
bet5_vector[i] <- win5
bet20_vector[i] <- win20
bet50_vector[i] <- win50
}
t5 <- table(bet5_vector)
t20 <- table(bet20_vector)
t50 <- table(bet50_vector)
print(sprintf("Apuesta = 5 -> E = %d, F = %d, Prob = %f", t5["1"], t5["0"], t5["1"]/(t5["1"] + t5["0"])))
## [1] "Apuesta = 5 -> E = 54, F = 4946, Prob = 0.010800"
hist(bet5_vector, main = sprintf("Éxitos y fracasos para apuesta = 5"))
print(sprintf("Apuesta = 20 -> E = %d, F = %d, Prob = %f", t20["1"], t20["0"], t20["1"]/(t20["1"] + t20["0"])))
## [1] "Apuesta = 20 -> E = 1217, F = 3783, Prob = 0.243400"
hist(bet20_vector, main = sprintf("Éxitos y fracasos para apuesta = 20"))
print(sprintf("Apuesta = 50 -> E = %d, F = %d, Prob = %f", t50["1"], t50["0"], t50["1"]/(t50["1"] + t50["0"])))
## [1] "Apuesta = 50 -> E = 1959, F = 3041, Prob = 0.391800"
hist(bet50_vector, main = sprintf("Éxitos y fracasos para apuesta = 50"))
Puede verse que la probabilidad de éxito va aumentando a medida que el monto de apuesta aumenta, ya que es más fácil ganar al ser menor la cantidad de juegos necesarios para pasar la barrera de los 200.
A work by CC6104